Phase diagram and momentum distribution of an interacting Bose gas in a 

bichromatic lattice 
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We determine the phase diagram and the momentum distribution for a one-dimensional Bose 
gas with repulsive short range interactions in the presence of a two-color lattice potential, with 
incommensurate ratio among the respective wave lengths, by using a combined numerical (DMRG) 
and analytical (bosonization) analysis. The system displays a delocalized (superfluid) phase at small 
values of the intensity of the secondary lattice V% and a localized (Bose glass-like) phase at larger- 
intensity Vi- We analyze the localization transition as a function of the height V2 beyond the known 
limits of free and hard-core bosons. We find that weak repulsive interactions unfavor the localized 
phase i. e. they increase the critical value of V2 at which localization occurs. In the case of integer 
filling of the primary lattice, the phase diagram at fixed density displays, in addition to a transition 
from a superfluid to a Bose glass phase, a transition to a Mott-insulating state for not too large V2 
and large repulsion. We also analyze the emergence of a Bose-glass phase by looking at the evolution 
of the Mott-insulator lobes when increasing V2. The Mott lobes shrink and disappear above a critical 
value of V2. Finally, we characterize the superfluid phase by the momentum distribution, and show 
that it displays a power-law decay at small momenta typical of Luttinger liquids, with an exponent 
depending on the combined effect of the interactions and of the secondary lattice. In addition, we 
observe two side peaks which are due to the diffraction of the Bose gas by the second lattice. This 
latter feature could be observed in current experiments as characteristics of pseudo-random Bose 
systems. 
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I. INTRODUCTION 

The interplay between disorder and interactions has 
been a long-standing challenge for condensed matter the- 
ory. In the absence of interactions a random potential 
can induce Anderson localization [lj, i.e. make all the 
single-particle eigenstates localized. In absence of disor- 
der, bosons on a lattice with repulsive interactions dis- 
play, for commensurate filling, a superfluid (SF) to Mott 
insulator (MI) transition as the repulsion is increased [H, 
with the superfluid phase displaying large density fluctu- 
ations and a gapless excitation spectrum, while the Mott 
phase is incompressible and has a gap in the excitation 
spectrum. If one considers both repulsive interactions 
and disorder, these two effects will compete: while disor- 
der makes the bosons localized, short-range repulsive in- 
teraction energy increases as the square of boson density 
and hence the total energy of the system is minimized by 
depleting the localized condensate towards a more uni- 
form density distribution. As a result, in a lattice Bose 
gas with short-range interactions a novel Bose-glass (BG) 
phase, non superfluid yet compressible, emerges between 
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the superfluid and the Mott-insulator 0. From the exper- 
imental point of view, it is possible to realize a system of 
bosons in a random potential by placing 4 He in porous 
media such as Vycor, aeorgels or xerogels[3|,|4j], or by us- 
ing artificially disordered Josephson junction networks 
Experiments in porous media revealed that the critical 
exponents of the normal-superfluid transition in Helium 
were different from those in pure helium in the case of 
aerogels and xerogels. However, the aerogel and xero- 
gel structures can hardly be described by a short range 
correlated random potential. In the case of Josephson 
junctions, localization of vortices was observed, but be- 
cause of dissipation, this system cannot be treated as 
fully coherent. The phase diagram of disordered boson 
system has also been intensively studied by Quantum 
Monte Carlo simulations [|| 0, Mdm- While some conjec- 
tures made in Ref . [2| could be confirmed, it appeared that 
very large system sizes were required to obtain reliable 
results. Due to the theoretical difficulty of the problem, 
one approach has been to reduce the spatial dimensional- 
ity. In one dimension, it is known that in the absence of 
interactions all states are localized as soon as the random 
potential is non-zero [l(| [HI- Moreover, powerful spe- 
cific techniques are available to handle the interactions; 
this is the case e.g. of the bosonization technique [l2| or 
of the Density Matrix Renormalization Group (DMRG) 
method [13, [lj]. For the specific case of a one-dimensional 



Bose gas subjected to an uncorrelated disorder (in ab- 
sence of a lattice), the phase diagram has been obtained 
by Giamarchi and Schulz [l5[ , showing that while for zero 
interactions the system is always localized, for nonzero 
values of the repulsive interactions a superfluid phase is 
possible at small values of disorder. Ref.[l5[ also pre- 
dicted that the non-superfluid (Bose-glass) phase of an 
interacting Bose gas is expected to differ markedly from 
the non-interacting Anderson-localized (AG) phase, e.g. 
the density profile of a Bose glass phase is rather uniform, 
in contrast to the highly inhomogeneous density profile 
of an ideal Bose gas in a disordered potential where all 
the particles occupy the lowest single-particle localized 
orbital. The phase diagram of a disordered, interacting 
Bose gas in one-dimension has been the subject of several 
numerical investigations by quantum Monte Carlo meth- 
ods strong coupling expansions fill, a nd Density- 
Matrix renormalization group approaches [l9| , that have 
established the existence of a Mott insulating phase sep- 
arated from the superfluid phase by a Bose glass phase 
for disorder not excessively strong. For stronger disor- 
der, these numerical studies have established that only 
the Bose glass and the superfluid is present. Also, the 
existence of a superfluid dome in the phase diagram has 
been obtained for the incommensurate casefl9| . 

With the development of atom cooling and trapping 
techniques, studying the Mott transition of bosons has 
become experimentally feasible [201 ]. Moreover, recent 
experiments with ultracold atomic gases have realized 
a pseudo-disordered potential by superimposing two op- 
tical lattices with incommensurate ratio between their 
spatial periodicities [2l| in a regime where interactions 
are important. Experimentally it is possible to charac- 
terize the system by measuring the excitation spectrum, 
the momentum distribution and higher-order (e.g. noise) 
correlations functions, as well as by looking at the equiva- 
lent of transport behavior through the study of the damp- 
ing of large-amplitude dipole oscillations [22| . 

While the experiments performed with a bichromatic 
lattice were focused on a regime where the lattice acts 
as a disorder potential, the physics of a bichromatic 
lattice is much richer, and the aim of this work is to 
describe the different possible phases of an interacting 
Bose gas subjected to such lattices. In the absence of 
interaction, the Schroedinger equation in a bichromatic 
potential treated in the tight binding approximation is 
known as the Harper model or the "almost Mathieu 
problem" and has been extensively studied by solid state 
phvsicistsf23l [2^| and mathematical physicists [25]. It is 
known to display a delocalized regime for weak incom- 
mensurate potential, and a localized regime for strong 
incommensurate potential, the two regimes being related 
by a duality transformation. In the limit of infinitely 
strong repulsions among the bosons (the so-called Tonks- 
Girardeau regime), the problem can be solved by map- 
ping to an ideal spinless Fermi gas subjected to the same 
external potential [261 ]. In particular, the model displays 
the same localization-delocalization threshold as in the 



incommensurate case 



AG 


BG 




7 


SF 





U/t 



commensurate case 









AG 




BG / 






• / 
/ MI 










SF 





U/t 



FIG. 1: Schematic representation of the expected phase di- 
agram for a Bose gas subjected to a bichromatic potential. 
"AG" is the Anderson-localized inhomogeneous phase, "BG" 
is the Bose glass phase, "SF" is the non-localized superfluid- 
like phase (ie displaying power-law decay of the phase-phase 
correlation function), and "MI" is the Mott insulator phase. 
The "?" sign stands for the region which need to be numeri- 
cally investigated. 



nonintcracting case. However, the momentum distribu- 
tion of the Tonks-Girardeau bosons is not directly related 
to the one of the spinless fermions, and for the specific 
case of the bichromatic potential it has been studied in 
[27l ]. The case of spinless fermions (or hard core bosons) 
with nearest neighbor repulsion was studied in [28[ . We 
focus here on the regime of intermediate repulsive in- 
teraction strengths. In the case of commensurate filling 
of the primary lattice and for A = a Mott-insulator 
phase is expected to occur at large values of interac- 
tion strengths U c /t ~ 3.3 [30| . In the disordered case, 
this Mott insulating phase competes with the localized 
phase, and is expected to induce a Bose glass intermedi- 
ate phase. 

A Bose gas subjected to a quasiperiodic potential with 
of finite interaction strengths has been previously stud- 
ied by Roth et al. [3l[ by exact diagonalization on a very 
small system and by Roscilde [32j in the case of a spe- 
cific choice of the height of the secondary lattice. In the 
present article, we use a combination of density matrix 
renormalization group (DMRG) methods and low-energy 
bosonization techniques to infer the phase diagram of the 
gas at varying height of the secondary lattice and interac- 
tion strengths, both for the case of integer and noninteger 
filling of the main lattice. The schematic summary of the 
known limits of the phase diagram is presented in FigfT] 
One of our aim is to see how the Mott lobes are modified 
by the presence of the secondary lattice in the commen- 
surate case and to establish a phase diagram for both the 
commensurate and incommensurate case. We also com- 
pute the momentum distribution of the gas, which is one 
of observables experimentally accessible. 

The paper is organized as follows: In Sec. II we intro- 
duce the model and the respective physical observables 
and give the low-energy description of the system via 
bosonization approach. Sec. Ill describes the numerical 
DMRG method. The results for the phase diagram both 
for non-integer and integer filling at varying the strength 
of the second lattice are given in Sec. IV. Here also the 
evolution of the Mott-lobes with pseudo-disorder is given. 
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In Sec.V we analyze the momentum distribution function 
and describe its characteristics for a weakly interacting 
Bose gas within perturbation theory in the strength of the 
second lattice. In Sec.V the dependence of the Luttinger 
exponent on pseudo-disorder is also determined. Finally, 
in Sec. VI we give a summary and the conclusions. 

II. MODEL 



where L is the length of the chain and E the ground state 
energy. 

We also evaluate the momentum distribution as the 
Fourier transform of the one-body density matrix, 

n(q)^Mj2^ m)a (b}b m )- (5) 

I in 

with a = n/ki being the primary lattice spacing and Af 
a normalization constant. 



We consider a one-dimensional Bose gas at zero tem- 
perature subjected to a bichromatic lattice potential 
V(x) = Vx sin 2 (fcix) + V 2 sm 2 (k 2 x): 



dx^\(x) (~^-V 2 + V(x)) i, b {x) 
\ 2m J 



+ o / dx^l{x)tpl{x)ip b (x)^ b (x), (1) 



H 



where ipb(x) is the bosonic field operator, m is the atomic 
mass and g represents the contact interaction. In the 
case where the main lattice is quite large, i.e. V\ > Er, 
where En = h 2 k 2 /2m is the recoil e nerg y, we can map 
the system on a Bose- Hubbard model 34] : 



N 3ite 



2=1 

NsUes N aU es 

1=1 1=1 



(2) 



where 6], hi are bosonic field operators on the site i, t 
is the hopping amplitude, U is the on-site interaction, p 
is the chemical potential, N S i tes is the total number of 
lattice sites; the parameters U, t are related to those of 
the continuum model (P)(e.g. see Refs. [H, HI]). The 
effect of the second lattice is to induce a modulation of 
the on-site energies according to A,; = Acos(27tcm), with 
A cx V 2 is the relative str eng th of the second lattice, 
and the value of a = k 2 /k\ [35| has been chosen as a = 
830/1076 ss 0.77, being the same of the experiment in 
Florence [U. 

In order to characterize the different phases of the sys- 
tem, we evaluate the following observables: (i) the super- 
fluid fraction, 



N 2 

sites 

Ntir 2 



.^antiPBC 



EpBc) 



(3) 



where N is the particle number, and E, 4 n dd ^ is the 
ground state energy with (anti)periodic boundary con- 
ditions, and (ii) the compressibility, \ — (l/L)dN/dfi, 
i.e. 



A. Low-energy properties and bosonization 

We focus now on the regime A <C 2t, which is expected 
to be nonlocalized [H, HE], [36[. In these conditions, we 
describe the one-dimensional interacting bosonic fluid as 
a Luttinger liquid, using a low-energy hydrodynamic de- 
scription [13, [33, [38| . In particular, the system is char- 
acterized by a slow, power-law decay of the phase-phase 
correlation function (hence the denomination of "super- 
fluid phase" ) with an exponent that depends on the in- 
teraction parameters. 

The low-energy Hamiltonian for the fluid can be writ- 
ten as 111, [33| : 



Ha = hj d 4 (V ^ (l))2 +v s K{^Ii{x)f 



(6) 



This Hamiltonian is a standard sound wave one in which 
the fluctuations of the phase 4>{x) represent the phonon 
modes of the density wave as given by 



p(x) = \po--V(f>(x)] V 

7T » — ' 

p= — 00 



o i2p(irp x-<j>{x)) 



(7) 



where po is the average density of particles. The 
field 6(x) — 7r J dx'H(x'), is conjugate of (f>(x), 
[~V(f)(x), 6{x')] = i5(x — x') and represents the phase 
of the supcrfluid. The parameters K and v s used in (J5|) 
are related to the microscopic compressibility and super- 
fluid density through the relations Kv s = irp s /m and 
v s / K — I/ttx- In the case of contact interaction between 
bosons g5(x) and in absence of the lattice potential, the 
Luttinger parameters v s and K are obtained by the ex- 
act solution of the Lieb-Liniger model [Hi: v R K = 
as follows from galilean invariance 



and = £ in the 

K 7T 



the 



weak coupling limit, while — ( 1 
strong coupling, po — N/L being the particle density. 
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The Hamiltonian ^ is an effective low-energy theory [J 
and provided that the correct values of the parameters 
v s ,K are used, all long wavelength properties of the cor- 
relation functions of the system then can be obtained 
exactly. In the g = 00 limit, i.e. for hard-core bosons 
one obtains K — 1 as for free spinless fermions while the 
free bosons limit would correspond to K — > 00. 

In the low-energy hydrodynamic description the 
bosonic field operator can be represented as 



L[E{N + 1) + E(N- 1) - 2E(N)}, (4) 



ip B {x) 



Vp( 3 



(8) 
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The corresponding one-body density matrix G(x, x',0) — 
(^(^OV'sCa-')) in the long- wavelength limit can be 
computed [L2[ and has a power-law decay given by ~ 
l/\x — x'\vk in the limit of the system size L — > oo. 
Notice that the knowledge of the compressibility and of 
the one-body density matrix offer two independent ways 
of extracting the Luttinger exponent 29] K . 



B. Perturbative treatment of the quasiperiodic 
potential 

For the model ([6]) we are interested in the effect of a 
bichromatic lattice potential V(x) = J2i=i Vi cos(2kix). 
We will work in the limit where the strength of both 
potentials are small with respect to the bandwidth, so 
that bosonization is applicable. Then, each component 
Vi of the potential couples to the density and adds a term 
to the Hamiltonian (J6j) which reads: 

Hu — Vi J dx cos(2kix) p(x) 

= / dxcos[{2-KppQ ±2ki)x-2p<t>{x)$) 

P— — 00 

Since the field <f>{x) is a slowly varying function on the 
scale of the interparticle distance, if oscillating terms re- 
main in the integral, they will average out, leading to a 
negligible contribution. Therefore, the Luttinger liquid 
(superfluid) behavior will persist provided that the filling 
is not commensurate i.e. none of the two the commensu- 
rability conditions ppo ± hi/ir € Z are satisfied. 

For commensurate fillings i. e. when one of the two 
commensurability conditions is met, the periodic poten- 
tial changes the simple quadratic Hamiltonian ([6]) of the 
Luttinger liquid into a sine- Gordon Hamiltonian which 
describes the Mott transition as a function of interaction 
strength [121 ]. Indeed, under the renormalization-group 
(RG) flow, the operator |9]) is irrelevant for K > K c = 
2/p 2 and relevant for K < K c , thus implying a Mott- 
insulator phase at K < K c . As K is decreasing when 
interactions are made more repulsive, this means that 
the Mott state is obtained when repulsion exceeds a crit- 
ical value U c . In the case where the Mott insulator is 
obtained for p — 1, in the regime of K < K c , none of the 
terms associated with the second potential (which is in- 
commensurate) can become relevant. Therefore, in that 
case, for K > K c , the Luttinger liquid is stable, and no 
Bose Glass phase can be created by the other potential 
in the vicinity of the Mott Insulator superfluid transition 
in the regime where bosonization is applicable. This jus- 
tifies the shape of the phase diagram of Fig. [1] for the 
commensurate case. The renomalization group analysis 
shows that the transition from the Mott insulator to the 
superfluid belongs to the Kosterlitz-Thouless universality 
class fl2l]. Note that the term (|9]) has been derived here 
for a weak lattice potential, but it appears also in the 



opposite limit of a strong lattice potential if the filling is 
commensurate, showing that the two limits are smoothly 
connected fl2^. 

A different situation occurs in the case of random dis- 
tributed disorder. As shown in Refs.[l5( the potential 
becomes relevant below the critical value K c — 3/2. Be- 
low such value the system lies in a Bose-glass phase with 
exponentially decaying Green's function on the scale of 
the localization length. A detailed RG analysis for the 
case of a generic quasiperiodic potential was given in 
Refs.[4Tl.l42l]. There it was shown that in the case where 
the quasiperiodic potential has a nontrivial, dense Fourier 
spectrum, the critical value of K c can be actually smaller 
than the value K c = 2, the deviation from K c = 2 be- 
ing related to the distance of 2irpo to a harmonic of the 
Fourier transform of the potential, thus interpolating be- 
tween the two-color potential and the fully random case. 

If we now consider the phase transition between the 
Mott state and the superfluid, not as a function of in- 
teraction, but as a function of particle density or as a 
function of the chemical potential, it is well known that 
in the absence of the secondary lattice potential, this is a 
commensurate-incommensurate (C-IC) transition [l2l l33l . 
[39l . |40| . At the transition, the scaling dimension of the 
operator cos 2(f) associated with the main lattice poten- 
tial must be 1, which yields K c = 1. Turning on a sec- 
ond, weak lattice potential incommensurate to the first, 
we see that the problem is reduced to free fermions in a 
bichromatic lattice. The rigorous results on the Harper 
model[25j then indicate that for a potential which is small 
compared with the bandwidth, the states are not lo- 
calized by the incommensurate potential. Therefore, a 
weak incommensurate potential cannot turn the super- 
fluid state formed by doping the Mott insulator in a Bose 
glass state. Again, this is at variance with the effect of 
the random potential, which would immediately localize 
the particles as the Mott gap closes. With model @, in 
the limit of very strong repulsion U 3> t, and for a filling 
slightly below one particle per site, we can also use the 
Harper model mapping to predict that the Bose glass to 
superfluid transition will happen when A = 2t. Thus, in 
the phase diagram at fixed U, and varying t, we expect 
that wings of Bose glass phase will be obtained for suf- 
ficiently small t. Summarizing the results for the Mott 
transition as a function of chemical potential and inter- 
action, we expect in the two-color potential a scenario 
similar to the scenario 2(c) in 2], i. e. that near the tip 
of the Mott lobe, there is no Bose glass phase in the case 
of the two-color potential, provided that the incommen- 
surate potential is small compared to the bandwidth. 



III. NUMERICAL METHOD 

In order to determine the ground state properties of 
the interacting Bose gas in the bichromatic lattice, we 
use the Density Matrix Renormalization Group (DMRG) 
method [II [l4j. The DMRG is a quasi-exact numerical 
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technique widely employed for studying strongly corre- 
lated systems in low dimensions. Based on the renor- 
malization, it finds efficiently the ground state of a rela- 
tively large system with quite high precision. Recently, 
the DMRG has already been used to study the ID disor- 
dered Bose-Hubbard model[43[. 

We consider a system with periodic boundary condi- 
tions and use first the infinite-size algorithm to build 
the Hamiltonian up to the length L, then we resort to 
the finite-size algorithm to increase the precision within 
many sweeps. In principle the Hilbert space of bosons 
is infinite; to keep a finite Hilbert space in the calcu- 
lation, we choose the maximal number of boson states 
approximately of the order 5(n), varying n max between 
n max = 6 and n max = 15, except close to the Anderson 
localization phase where we choose the maximal boson 
states n max = N. The number of eigenstates of the re- 
duced density matrix are chosen in the range 80 — 200. To 
check the error produced by truncating the boson space, 
we have repeated the calculations at varying n max in the 
range 5(n) and 10(n), without observing substantial dif- 
ference in the ground state energy. To test the accuracy 
of our DMRG method, in the case U = or for finite U 
and small chain, we have compared the DMRG numerical 
results with the exact solution obtained by direct diag- 
onalization. For larger system (N S i tes > 10), we have 
checked the convergence of the ground state energy by 
varying the number of truncated eigenstates, estimating 
that in the region of the superfluid-Mott insulating phase 
the errors are of the order 10~ 6 . The good convergence 
of the algorithm is also tested by the coherence of the 
results obtained from different observables as the Mott- 
insulator density plateaus and correlation functions. 

The calculations are performed in the canonical ensem- 
ble, i.e. at fixed number of particles N. The chemical 
potential is determined by the evaluation of the energy 
required to add or subtract a particle to the ground state, 
i.e. n p = E(N+l)-E(N) and fi h = E(N)-E(N-l)\v^. 
In this way we may obtain the phase diagram in the grand 
canonical ensemble. In order to find the superfluid den- 
sity and the compressibility at varying chemical poten- 
tial, we performed several calculations at varying particle 
numbers. For the determination of the phase diagram we 
have chosen N S i tes = 20, while the correlation functions 
have been calculated using a larger chain N sit es = 50. 



IV. PHASE DIAGRAMS 

We have determined the phase diagram in two situa- 
tions. First, we have analyzed the effect of interactions 
on the localization/delocalization threshold with respect 
to its noninteracting value A = It obtained from the 
Harper model [23| or equivalently for the hard-core Bose 
gas. Secondly, we have analyzed the effect of disorder on 
the Mott-insulator lobes Q. 



A. Localization/delocalization transition 

1. Incommensurate filling: case (n) = 1/2 

By the calculation in the canonical ensemble of the su- 
perfluid fraction and of the compressibility, we have eval- 
uated the phase diagram in the plane (A/t,U/t). This is 
illustrated in Fig. [5] (upper panel) by showing the contour 
plot of the superfluid fraction obtained for N S it es = 20. 
In the case of non-integer filling only two phases are 
present: a superfluid phase (f s ^ 0) at small values of 
the secondary lattice height A (bottom-left), and a Bose 
glass phase (f s = 0) at large values of A for A > U (top- 
left). At U = the transition occurs at the expected 
critical value A/t = 2. We see that at intermediate val- 
ues of the interaction strengths U the critical value of 
A c /t increases and there the superfluid region extends 
in a large dome. A similar behavior is observed for a 
disordered Bose gas (l5| . 

2. Commensurate filling: case (n) = 1 

The phase diagram for the integer filling is given in 
Fig. [5] (lower panel) where are reported the superfluid 
fraction f s (main figure) and the compressibility gap 
(fi p - fi h )/t (inset) obtained for N mtes = 20. The Mott- 
phase which is characterized by a large compressibility 
gap emerges at the bottom right corner above the criti- 
cal value U c /t = 3.3 ± 0.2 for A = in agreement with 
Ref . [29l |30| . We observe that U c increases at increas- 
ing A, meaning that disorder energetically reduces the 
compressibility gap in the localized regime (see also the 
upper panel of FigJ3] below). A Bose-glass phase instead 
occurs in the region of the phase diagram A > [/(top- 
left). At U = the transition occurs at the expected 
value A c /t = 2. The critical value of A c increases with 
U at small U indicating a derealization by interactions, 
similarly to the true-disorder case. Finally, a superfluid 
phase emerges in the small U and small A region of the 
phase diagram (bottom- left). In our simulations it ex- 
tends in a large dome at intermediate U and A. The be- 
havior of the superfluid fraction and compressibility gap 
for small A and intermediate U seems to indicate a di- 
rect transition from the superfluid to the Mott-insulating 
state without passing into a Bose glass. Such conclusion 
seems physically reasonable if one takes into account that 
the bichromatic lattice potential acts as a quasi-disorder, 
i.e. is less relevant than true disorder. Anyway such 
conclusion should be supported by further numerical in- 
vestigation and finite size scaling of the compressibility 
and superfluid fraction. 

B. Mott-insulator lobes 

We have performed the calculation of the Mott- 
insulator lobes in the grand canonical ensemble. This is 
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FIG. 2: DMRG phase diagram for an interacting Bose gas 
in a two-color lattice, in the plane (A/t , U/t). Upper 
panel: superfluid fraction in the case of non-integer filling 
v = N/N sites = 0.5, with N = 10, iV s ; tes = 20. Lower 
panel: the supernuid fraction f a (main figure) and compress- 
ibility gap (n P — fir)/t (inset) in the case of integer filling with 

N = Nsites = 20. 



obtained by the estimation of n p and /i h for several values 
of particle numbers. At increasing strength of the second 
lattice we find that the Mott-insulator lobe with (n) = 1 
shrinks and finally tends to disappear for A ~ 0.5, as is 
illustrated in Fig. [3] (upper panel). In order to determine 
the Bose glass region we have also calculated the super- 
fluid density. Fig. [3J (lower panel) shows, for a specific 
choice of A, the regions of nonzero superfluid density as 
well as the regions of large compressibility gap (Mott- 
insulator phase) through the function f s + (fi p — fJ, h )/t. 
The intermediate (dark blue) region between the two cor- 
responds to the Bose glass phase. Notice that near the 
tip of the Mott lobe the superfluid fraction is nonzero, 
as illustrated in the inset of Fig. flower-panel), support- 
ing the direct superfluid to Mott-insulator transition, dis- 
cussed above. We also notice on Fig. [3] the presence of a 
Bose glass phase for t/U < A./2U, as expected from the 
strong coupling argument. 




FIG. 3: DMRG phase diagram for an interacting Bose gas in 
a two-color lattice, in the plane (n/t , t/U), for the first Mott 
lobe, for N = N S it e3 = 20. Upper panel: the shrinking of the 
Mott lobe at varying A/U=0 (solid line), 0.1 (circles), 0.2 
(squares), 0.5 (diamonds). Lower panel: contour plot of the 
function f s + (/i p - fi h )/t for A/U = 0.1. The inset shows the 
compressibility gap (fi p — fi h )/t and the superfluid fraction f B 
along the line jj,/U = 0.25. 




qa/ix 



FIG. 4: DMRG momentum distribution function in the su- 
perfluid phase at varying A/U (as indicated on the figure) 
and U — 2t, with N = N S it es = 50. Subdominant peaks are 
determined by the presence of the second lattice potential (see 
text). 



V. MOMENTUM DISTRIBUTION 

A. Side peaks of the momentum distribution 

The results for the momentum distribution are re- 
ported in FigdJ We note that besides the expected peak 
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of the momentum distribution at k — 0, there are other 
peaks at k = ±Q = ±2^(1 — a) related to the modulation 
of the on-site energy in Eq. ([2j) • The origin of these peaks 
can be understood by considering first non-interacting 
bosons. We will begin by discussing the continuum limit, 
and then the lattice case. If we approximate the irra- 
tional number a [35| by a rational number p/q, in the 
potential V(x), we can apply Bloch's theorem and write 
the boson annihilation operator as: 



4>b(x) 



(10) 



k (3=1 



where k is the quasi-momentum of the boson, and (3 is 
the band index. Bose condensation will then occur in 
the lowest quasi-momentum state of the lowest band (we 
chose j3 = 1 for this band). The functions Lpf\x) are 

periodic of period qa, i.e. (pf\x) = f^\x + qa). Using 
this property one finds that in the Bose Condensed state, 
(%(x + qa)ipb(x' + qa)) = (^l(x)ip b (x')) . As a result, 
after averaging over x, the function {ipl(x + r)ipt(x)) be- 
comes a periodic function of r of period qa. Using Fourier 
transformation, we conclude that the states of momen- 
tum (2n/a)(l/q ± m) present a macroscopic occupation 
number. If we turn to perturbation theory, in the limit 
of A <C t, we find that the perturbed wavefunction at the 
lowest order is given by 



(i),o 



(,;) 



(11) 



£ 

Q,m 



("0,0/ \ 

/* [X J (^ m) -°(^|U 2 cos(2afc 1 x)|^ 1) -°(,)), 



where i fi\ m ^ '°(x) are the solutions of a Mathieu 
equation [45j for the potential V\cos(2k\x) and Eq. n is 
the dispersion of the n— th band for momentum Q. The 
matrix elements of perturbation are non zero only when 
Q = Q± = (2n/a)(a ± m) (me Z). 

The momentum distribution is then given by n(p) = 

| J dxe wx vi=n [x)\ 2 and using Eq. tfTTj) we find that it dis- 
plays two peaks: 



n(f) ~ WMt 



E 1 * 

s=± ^ 



Vo 



En 



\^(P + Qs)\ 2 . (12) 



1/4 



where)!!] pg(p) oc e - p2 / p ° and p = f 

In an analogous way we can proceed to derive the ex- 
pression for the momentum distribution on the lattice. 
The perturbed boson annihilation operator is then: 



b k = 



£ 

5=± 



A 



-2i(cos((/c + Qs)a) — cos(fca)) 



so that the largest occupation number will be found for 
k = 0, and again at k = Q± (modulo the reciprocal 
lattice vector). The physical interpretation of the extra 



peaks is therefore that the ground state wavefunction is 
diffracted by the quasiperiodic potential thus creating 
peaks at multiple harmonics of 2ira/a (modulo a vector 
of the reciprocal lattice). 

Let us now turn to the case of weakly interacting 
bosons. If the repulsion U is not too large, we can 
still begin by diagonalizing the non-interacting Hamil- 
tonian, and treat the interaction within Bogoliubov ap- 
proximation or numerically solve the Gross-Pitaevski 
equation . Since Bose condensation is obtained in the 
lowest band, it is reasonable to neglect the contribution 
from the higher bands. Moreover, since the states that 
are important for the low energy properties are those 
with quasi-momentum near zero, we can neglect the de- 
pendence of <p^\xi) on k. This gives us the following 
expression for the boson annihilation operator [47| : 



ip^\xi)bi 



(14) 



where bi = jjr/2^k elkXi bk,i- Injecting this approxima- 
tion in the full Hamiltonian, we obtain an interaction 
term which has the same period q as the potential A^. 
This gives rise to new umklapp processes, but since we 
are only interested in the states of momenta close to zero, 
we can neglect them. Then, the theory describing the b 
bosons becomes identical to the one describing bosons in 
the absence of incommensurate potential, albeit with a 
dispersion fixed by the band structure and an interaction 

u,f f . = uY&\$\xi)\ A /<i. 

The single particle density matrix is: 



(15) 



and thus the effect of the periodic potential is only seen in 
the appearance of the factor {fo (xi))*(pg(xj). Using 
the bosonization technique to compute the single particle 
density matrix we finally find that: 



(16) 



By Fourier transforming the above expression, we recover 
power law peaks in the momentum distribution with ex- 
ponent [1/2K — 1] for all the wavevectors that are mul- 
tiples of 2ir/qa . Based on the previous perturbation 
theory, we expect that the two subleading peaks will be 
found at k = (2ir/a)(m±p/q). Moreover, the exponent 
should be identical to the one found for q — 0. We also 
remark that if the peaks were produced by the terms 
e i2trp x e i(8-2<p) m £ ne ex p ans i on f the boson annihila- 
tion operator {§}, their position would depend on the 
number of particle per site, and their height would be 
independent of the strength of the incommensurate po- 
tential. Moreover, these terms give in real space a corre- 
lation function of the form (\x - x'\/a)- {2K+1/2K) with 
an exponent that is always larger than two. As a result, 
the Fourier transform of this term would not diverge as 
k —* (27r/a)(m ±p/q), than a cusp could be obtained. 



FIG. 5: Fourier transform of DMRG momentum distribution 
function in the superfluid phase with A = 0.5C7, U — 2t and 
N = N S ites = 50. The main peak and the subdominant one 
decay with a power-law exponent consistent with [1/2K— 1] ~ 
0.85 for q sufficiently close to and 27r(l — a)/a as shown in 
log scale in the insets. 



We have checked that the height of the secondary peak 
increases quadratically with the strength of the incom- 
mensurate potential as expected from (JX3J) , that its posi- 
tion does not change with the filling, and that it possesses 
the same power law divergence as the peak obtained at 
k = 0. This is shown in FigJS] where the Fourier trans- 
form of the momentum distribution is displayed together 
with the power-law decay of the peak at q = and of the 
satellite peak in a log scale. 



B. Determination of the Luttinger exponent 



According to Eq. (fl6|) , in the superfluid phase the one- 
body density matrix pi(i,j) = (b\bj) can be used to ex- 
tract the Luttinger exponent K. This is particularly in- 
teresting because, even though bosonization techniques 
do not directly access to the localized phase, the fact that 
the Luttinger exponent K depends on the strength A of 
the pseudo-disorder indicates a first disruption of the su- 
perfluid phase towards localization. In order to analyze 
the DMRG data for the one-body density matrix, we take 
into account both the density modulation induced by the 
second lattice (entering explicitly in Eq. I]16p through the 
factors <po(xi)), and the fact that the calculations are 
performed on a system of finite length L. For the latter 
case, we use the results of the continuum model obtained 
by using the conformal field theory [38| for a system of 
length L and periodic boundary conditions. In essence, 
we fit the DMRG results by the following expression: 



PiO',0) = ^/n + S cos(2tt(1 - a)j + <f) ) 



d(ja\L) 



where no, 5 and <po are constants, K is the Luttinger 
parameter and d is the conformal length d{x\L) = 
— |sin(^)|. The results are shown in FigO By the fit 
we obtain that the Luttinger exponent K decreases at 



FIG. 6: One-body density matrix from DMRG data (A/17 = 
circles, A/U = 0.1 squares, A/U = 0.5 triangles) and from 
fit to Eq. (|17[) (dashed lines). The parameters used are U = 2t 
and TV = N szte s = 50. 



A/U 


K 


0. 
0.1 
0.5 


3.44± 0.03 
3.43± 0.04 
3.35± 0.06 



TABLE I: Values of the Luttinger exponent from the fit of 
the DMRG data to Eq.(TJ7Jl with the parameters of Fig[5] 
The corresponding \ IS of the order of 5 x 10 . 



increasing A, in agreement with the intuition that dis- 
order drives the system towards a more correlated, less 
superfluid phase. The corresponding values are reported 
in Table [J 

Another independent way to extract K is based on 
the determination of the ground state energy and com- 
pressibility x given by Eq.(j4]), by the relation K — 
hm\J psx/ m - We have verified that the values of K ex- 
tracted in this way are consistent with those of Table U 



VI. SUMMARY AND CONCLUDING 
REMARKS 

We have analyzed the phase diagram of an interacting 
one-dimensional Bose gas in the presence of a pseudo- 
disorder generated by a bichromatic lattice potential. 
Starting from a Bose-Hubbard model we have considered 
both commensurate and incommensurate fillings and we 
have found a rich phase diagram including, in addition 
to the superfluid and Mott phases, a Bose glass phase, 
localized but compressible. In agreement with the lim- 
iting cases of free and hard-core bosons described by an 
almost Mathieu problem, the transition towards the Bose 
glass phase is found at A/i > 2, the critical value of A 
being higher for bosons with finite interaction strength. 
This non-monotonic dependence of the critical height of 
the second lattice on the interaction strength could be 
observed in the experiments. We have also analyzed the 
shrinking of the Mott-lobes as a function of A and the 
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emergence of a Bose-glass phase in the (fi/U, t/U) plane. 
Finally we have characterized the superfluid phase by a 
static observable, the momentum distribution function. 
We have shown that satellites peaks emerge when the 
pseudo-disorder is not too strong and their interpretation 
within perturbation theory offer a good qualitative un- 
derstanding of their behavior as a function of the height 
of the second lattice. The central peak of the momentum 
distribution allows to determine the Luttinger exponent 
K, whose knowledge is useful to make predictions for 
further physical quantities. 

While the momentum distribution and the behavior of 
the side peaks could characterize the evolution of the sys- 



tem towards a Bose glass, a direct probe of a Bose glass 
phase and its distinction from a Mott-insulator could be 
provided by study of noise correlations or collective exci- 
tations. This will be left for future study. 
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